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Abstract. We examine the existence of shock profiles for a hyperbolic-elliptic system aris- 
ing in radiation hydrodynamics. The algebraic-differential system for the wave profile 
is reduced to a standard two-dimensional form that is analyzed in details showing the 
existence of heteroclinic connection between the two singular points of the system for 
any distance between the corresponding asymptotic states of the original model. De- 
pending on the location of these asymptotic states, the profile can be either continuous 
or possesses at most one point of discontinuity. Moreover, a sharp threshold relative 
to presence of an internal absolute maximum in the temperature profile -also called 
ZeI'dovich spike- is rigourously derived. 
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1. Introduction 

Fluid dynamics equations support a class of special and significant solutions known as 
shock waves. These describe physical phenomena consisting in an abrupt transition from 
one state to another, that, in the simplest setting, are described by a single jump connect- 
ing two different states and propagating, in first approximation, with constant speed. In 
formulas, these corresponds to special solutions of the system under consideration having 
the form 

w{x,t) = W{x - ct) W{±oo) = W±, 

where w denotes the state variable, W the shock profile, c the velocity and W± the 
asymptotic states. The difference [W] := W+ — W- and its length are generally used to 
distinguish between small, |W+ — W-\ — )• 0, and large shocks, \W+ — W-\ — )• oo. 

In the relevant case of the Euler equations for gasdynamics in one space dimension, 
obtained considering the conservation of mass, momentum and energy and neglecting 
higher order effect as viscosity, capillarity, and thermal conductivity, traveling wave solu- 
tions propagate with a velocity dictated by the Rankine-Hugoniot condition. Moreover, 
they have piecewise constant profiles with a single jump point, where physical entropy in- 
creases along the discontinuity. The persistence of the presence of waves when considering 
also higher order effects has been considered by several authors and in different context. 
Among others, let us quote the classical reference 0] where the existence of shock layers 
for a model taking into account both viscosity and heat conductivity. 
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The present article fits into the strand of radiation hydrodynamics (see |13lll6j ). Specif- 
ically, we consider the classic Euler equations coupled with an elliptic equation for an addi- 
tional variable describing the intensity of radiation; precisely, we deal with the hyperbolic- 
elliptic system 

' dtp + dx{pu) = 0, 

dt{pu) + dx{pu^ +p) = 
\dt{pE) + d^{pEu + pu-a~^d:,n) = Q, 
^dln = T\n-g{e)) 

where the variables p,u,p,E and 6 represent density, velocity, pressure, specific total en- 
ergy and temperature of the fluid under consideration. The additional variable n describes 
the specific intensity of radiation and the coupling is governed by the constants r, as and 
by the function g = g{9). 

The model ([l]) derives from a hyperbolic-kinetic system where the radiation is described 
by an additional variable for the photons density satisfying a transport equation with 
interaction kernel given by the Stefan-Boltzmann law. System ([T]) is obtained in the non- 
relativistic limit (speed of light tending to +oo) and the variable n is an average of the 
photon density. The elliptic equation for the variable n emerges from the assumption that 
the variable n is linked with the temperature by the integral relation 

n{x,t) = I g{e{y,t))K{\x-y\)dy 

with an exponential kernel K{s) = ^re"'^^. Details on the derivation of asymptotic 
regimes as the one described by ([T| can be found in [H [5l [12] . 

As usual, the pressure is considered a function of density and temperature, p = p{p, 6) 
and the specific total energy E is given by 

E = ]^v? + e 

where the internal energy e, in general, may depend on p and 6, e = e{p,6). Precisely, 
we concentrate on the case of polytropic perfect gases; this amounts in considering the 
functions p, e, g given by the formulas 

(2) p = Rpe, e = ^^, gie) = ae'' 

where 7 > 1 is the adiabatic constant (7 = 5/3 for monoatomic gases, and 7 = 7/5 for 
diatomic gases), R > the perfect gas constant, a the Stefan-Boltzmann constant and 
a = 4. 

Existence of shock profiles in the context of radiation hydrodynamics have been con- 
sidered for many years. In [71 [THl I13| . a formal analysis is performed concentrating on 
specific ranges for the asymptotic states relative to the size of the shock and the strength 
of the radiation effect. More recently, rigorous proof has been presented: 
- for a simplified model, known as the Hamer model, consisting of a single scalar conser- 
vation law coupled with an elliptic equation, mimicking the form of the system ([T]) (small 
and smooth shocks [T5], general and possibly discontinuous shocks [8} \9\ flO]): 
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- for systems, both for the special form ([T]) (small and smooth shocks {Xlj, large and dis- 
continuous shocks [2]) and for general hyperbolic-elliptic systems (small shocks and linear 
coupling [S], small shocks and nonlinear coupling [TU]). 

Here, we contribute to the exploration of system ([T]) showing existence of shock profiles 
for all of the possible regimes: small, medium and large size. Apart for the requests on the 
form of constutive functions p, e and g, the unique assumption is that the asymptotic states 
connected by the shock profile determine an admissible shock wave for the corresponding 
reduced hyperbolic system. 

Theorem 1.1. Assume ^ with > 0, 7 > 1, a > 0. Let the triples {p±,u±,d±) and the 
constant c GM be such that the Rankine-Hugoniot conditions 

c[p] = [pu], c[pu] = [pU^ +pI c[pE] = [pEU + pU] 

are satisfied, together with the entropy condition for 1 -shocks 

C+ < C<U+<C+y/^9^ 

Then system ([T]) supports a shock wave, that is solution of the form 

{p, u, 6, n) = {p, u, 9, n){x — ct) 

such that {p,u,9,n,dxn){±oo) = {p±,u±,9±, g{9±),0). 

The profile is unique up to translations and possesses at most one discontinuity point; 
if this is the case, all of the variables p, u and 9 have a jump at such point. 

Once the existence of the shock waves is estabilished, it is meaningful to determine 
qualitative properties of the profile and, precisely, its smoothness and its monotonicity. 
To this aim, we introduce the velocities U± := u± — c, giving the speed of propagation of 
the asymptotic states relatively to a reference frame joint with the shock. The following 
statement concerns with the regularity of the profile, giving sufficient conditions for the 
presence of a jump point, and with the monotonicity of the component of the profile. Both 
they are stated in term of the key parameter given by the ratio between the velocity jump 
of the shock, i.e. — U-\ = — U-\, and the average velocity given by the arithmetic 
mean (U^ + U-)/2. Precisely, we consider the parameter 

S := ^^'^ where Uc := ^ ^~ , 

so that the measure of the shock size is always normalized with respect to the average 
speed of propagation of the asymptotic states relatively to the shock velocity. 

Theorem 1.2. Assume ^ with i? > 0,7 > l,a > 0. 



i. If J < 2, the profile of the shock wave determined in Theorem 1.1 is discontinuous in 
the components p, u and 9 if 

2(7 - 1) _ 

7 
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ii. The component p and u of the profile of the shock determined in Theorem are 
always monotone, while the temperature profile is monotone if and only if ^ < 3 and 

7-1 
S < 5s,,ue := — . 

// this condition is not satisfied, the temperature profile has a single internal absolute 
maximum point that can he either attained at the jump point or at some point of regularity. 

As a consequence of the entropy conditions, the value of 5 belong always to the interval 
(0,2/7). For 7 G (1,3), there holds < (5spike < (^jump < 2/7, hence, depending on the 
range of values of the parameter (5, we can distinguish three different type of temperature 
profiles: monotone and continuous, non-monotone and continuous, non-monotone and 
discontinuous. 

The smallness of the parameter 5 is not equivalent to the one of the shock stength 
|(p_l_, C/+, ^+) — (p_, {/_, 6*-)!. Indeed, from the Rankine-Hugoniot relations, we infer 

and thus the small shock assumption corresponds to the requests 

M2<5«|[n]|i/2«i. 



As showed in [11^110). in this limiting regime, the profile is smooth. Theorem 1.2 gives the 
necessary condition 5 < (5jump for such regularity. An additional condition is also obtained 



during the proof, see (26), and not reported here because of its scarce readability. 



With respect to previous rigorous results on radiative shocks Theorems |1.1| and 1.2 



can be considered as conspicous improvements. By restricting the attention to a specific 
system, we can implement a strategy closely resembling the one considered in [U [TO] 
passing from a local study for small shocks to a global analysis for any possible shock size. 
With respect to [2\, where the case of large shocks is considered, in addition to the removal 
of the size assumption, we reduce drastically the assumption on the data of the problem. 

As usual, the proof is based on the analysis of the system of ordinary differential equa- 
tions obtained by looking for traveling wave solutions to ([T]). Taking advantage of the 
conserved quantity, such fifth-order system can be reduced to a second-order one. The 
main original ingredient in the proof resides in the choice of the variable selected for the 
reduction. The analysis in [21 [71 [8l [TI] is based on the study of the second-order differential 
equation solved by the velocity variable u; differently, in [9l[l0], it is preferred to consider 
the second-order equation for the radiation variable n. Here, we decide for the somewhat 
intermediate approach given by considering the first-order system for velocity u and radia- 
tion n. Such system is transformed in the standard form, that is common to the one valid 
for the Hamer model, and the subsequent analysis consists in a detailed study of such a 
prototypical system, for which the existence of a heteroclinic orbit, possibly discontinuous 
is proved (see Proposition 3.4). Theorems [1.1[ and 1.2 descend from Proposition 3.4 and 



the construction of the heteroclinic orbit for the reduced system. 
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The paper is organized as follows. In Section 2, we show how to reduce the original 
problem for the radiation hydrodynamics model ([T]) to a standard form. For completeness, 
we also present how the same strategy operates in the case of the Hamer model. Section 3 
is devoted to the analysis of the system in standard form and to the construction of hete- 
roclinic orbits connecting the two singular points. In Section 4, we resume to the original 
model, interpreting the construction of the previous Section in term of the hydrodynamics 
variables. Finally, we draw conclusions in Section 5, pointing at some possible appealing 
directions for future research in the subject. 

2. Scaling down to the reduced system 

Due to the conservation structure, the system of ordinary differential equations for the 
profile of the radiative shock under investigation can be reduced to a first order system 
for a two-dimensional variable. The aim of this Section is to show that such a reduced 
system has the form 

dx 1 G(x) — y dv 1 , o, 

for an appropriate function G and positive parameter i'. The same procedure can be 
applied to a simpler system given by the coupling of a scalar conservation law with an 
elliptic equation. For pedagogical reasons, we first present this simpler situation and 
postpone the discussion of the Euler system with radiation ([T]). Since the two derivations 
of ([3]) are independent, the reader may skip the first part and move directly to the case 
of radiation hydrodynamics. 

Hamer model. Given smooth functions / and g, let us consider the hyperbolic-elliptic 
system, known as Hamer model (see l6j and descendants) 

dtu + {f{u) - dxn) = 0, 
d'^n — n + g{u) = 0. 

For convex fluxes /, the hyperbolic scalar conservation law obtained by disregarding the 
coupling term dxTi possesses shock waves connecting states u± if and only if < m_. 
Such waves have a piecewise constant profile with a single jump point and propagates with 
a speed c that is given by the Rankine-Hugoniot condition, viz. 

^ ^ [m\ ^ f{u+)-f{u.) 

[u] U4. — li_ 

In analogy, we look for traveling fronts {u,n) = {u,n){x — ct) that solves Q and satisfy 
the asymptotic conditions 

u(±oo) = u± n(iboo) = g{u±), dxn{±oo) = 0. 

Apart for the requests ti+ < u_ and c given by ([s]), we also assume 

(6) [g] := g{u+) - g{u^) < 0. 

Setting ^ = X — ct, the couple {u, n) satisfies an algebraic-differential system 

(7) rn = f{u)- f{u±)-c{u-u±), ^ = m, =n- g{u), 



(4) 
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to be interpreted as a two-dimensional dynamical system along a surface in the three- 
dimensional space (u, n, m). 

There are different possible ways to deal with Q, corresponding to different choices 
of unknown. The approach in |15l |8] consists essentially in deriving a second order dif- 
ferential equation for the variable u, while the strategy used in O [10], dictated by the 
better regularity properties of the variable n, corresponds to consider the couple {n,m). 
Differently, we deal with ([7]) in a way that is somewhat in between the two and we look 
for a reduced system for the couple (n, n). This method avoids the necessity to invert the 
relation m = f{u) — f{u±) — c{u — u±), but it has the drawback of the possible presence 
of jumps of the trajectory. 

Differentiating the first equation, we get 

dm I J,,, \ du 

so that system Q is completely described by the reduced system 
du n — g{u) dn 



di f'{u)-c' di 



f{u) - f{u±) - c{u - u±). 



In the case of Burgers equation, i.e. f{u) = 2^ ^ c = (n+ + u-)/2, 
du n — q(u) dn 1 , 

Introducing the variable {x, y) related with (ti, n) by 

u = ^|[M]x + n_ +n+| and n =]^\[g\y + g{u-) + g{u+)^, 
where [u] = u+ — u^ < 0, system ([s]) takes the form ([s]) where 

2 [5] ^ ' [9] 

The function G and the parameter have to be regarded as functions of the asymptotic 
states u± or, equivalently, of the values [u] and c. By definition, for any choice of u± there 
hold G(±l) = ±1. 

Example 2.1. For a linear term g, viz. g[u) = au for some a > 0, the expressions for v 
and G are 

V = —— [u] and G{x) = x. 
In the quadratic case, g{u) = au^, a > 0, we obtain 

j/ = --^M and G{x)=x-]^^{l-x^). 
Aa c 4 c 

where c = + u-). In this case, both u and G are homogeneous with respect to the 

couple {[u], c) and thus system ([s]) does not vary if the ratio [u]/c is kept fixed. 

In general, for g{u) = au" for some exponent a > 0, the value [g] is homogeneous of 

degree a with respect to [u] and c, and, as a consequence, the parameter 1/ is a homogeneous 

function of degree 2— a with respect to the couple {[u], c) and G is homogeneous of degree 0. 

Explicit expressions for the corresponding i' and G does not seem particularly significant. 



Solutions to system ^ may possess eventual discontinuity. At such point, y is contin- 
uous (as a consequence of the continuity of n) and, denoting by and Xr the values of x 
at the left and at the right of the jump point, respectively, there hold 

X£ < Xr and Xi + Xr = 0. 

Solutions to Q with the above jump conditions will be analyzed in the Section [sj 

In the remaining part of the present Section, we show how system ([s]) emerges also 
when looking for traveling wave solutions fot the radiating Euler system ([T]). 

Radiation hydrodynamics. Next, we consider the hyperbolic-elliptic system ([T]) de- 
scribing the evolution of a compressible fluid under the effect of radiation. Setting = 0, 
the first three equations in system ([T]) reduce to a standard hyperbolic model for com- 
pressible fluids. Such system possesses shock solutions, i.e. traveling wave solutions 
W = W{x — ct) with piecewise constant profile W of the form 

for appropriate states W± and parameter c. The shock speed c is related with the states 
W± by means of the Rankine-Hugoniot condition. Precisely, setting U = u — c, there hold 
( P-U- = p+ U+, 
p- U^+p_= p+ ul+p+, 

P-U- (lu^+e-^ +p-U-=p+U+(^^Ul + e+^ +P+U+ 
In the case ([2]), if p±U± 7^ 0, the values 9± satisfy the linear system 

u^0_ _ = i u-u+[u], e--e+ = (c/+ + U-)[u\, 

hence, they can be written as explicit functions of U± 

where Uc ■= ^{U- + Note also that, if the jump [u] is strictly negative, then the 

jumps [p] and [6] are both strictly positive. 

Still taking advantage of ([2]), entropy conditions in the case of 1-shocks require U±,6± 
to be such that 

(9) <u-, < [/+ < 

which implies f7+ < C/_, since [9] and [u] have opposite sign. 
For later use, let us set 

(10) A:=p±U±, B:=p±Ul+p± C := p± t/± Q + e±^ + p± C/±. 
or, equivalently, 

A:=p±U±, B:=2±l^ ([/_ + [/+) C:= AU.U+. 

27 2(7 - 1) 

having used Q. 
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We are interested in determining the existence and the internal structure of travehng 
wave solutions 

u, 6, n) = {p, n, 9, n){x — ct) 

solutions to ([T]) such that {p,u,9,n,dxn){±oo) = {p±,u±,9±, g{9±),0). Plugging the 
traveling wave ansatz into the system, we get the system of differential equation 

^ipU) = 0, l(p[/2+p) =0, 

to be satisfied in the regions where the wave profile is smooth. 
The conservative form of the first two equations implies that 

pU = A and pU"^ + p = B 



all along the trajectory for any weak solution to (11) where A and B are given in (10). 
Thus, there hold 

(12) p = p([7) = ^, p = p{U) = B-AU 



and system (11) reduces to the algebraic-differential system 

(13) m = asf{U), ^=m, = r\n - {g o 9){U)) 



where, for A,B,C as in (10), 

(14) f{U):=-^AU^ + Ae{U) + BU-C, 

and e{U) = e{p{U), 9{U)) with p{U) and 9{U) are obtained from (U and ([12]). System ^ 
describes a two-dimensional dynamical system in the three-dimensional space {U,n,m) 
along the manifold S determined by the algebraic relation m = asf{U). Using the special 
form for pressure p and internal energy e given in ^ , the function / can be rewritten as 

(15) f{U) = -K{U-U+){U-U-) where k := J^\, A. 

2(7 - 1) 

If the profile has a jump point at some point, the function n and its first derivative are 
forced to be continuous at such point, since the last equation in (13) forces the second 
derivative of n to be bounded and measurable. As a consequence, denoting by p£,Ue,9£ 
and pr,Ur,9r the values at the left-/right-hand side of the discontinuity, the relation 
f{Ui) = f{Ur) holds together with the admissibility conditions that read as, in the case 
of 1-shocks, 

VRlOi <Ue, 0<Ur< y^RjTr 



Because of (15), jumps are admissible if and only if 

(16) Uc<Ue< 2Uc where Uc := ^ 

the corresponding value Ur being given by 2Uc — Ui. 
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There are different possible ways to deal with (13), corresponding to different choices 
of unknown. As for the Hamer model, two main approaches have been used: in [71 [TT| [2] 
a second order differential equation for the variable U is derived, in [SJ [TU] it has been 
considered the dynamics for the couple {n,m), because of the greater regularity of the 
variable n. Here, we proceed by deriving a system for the couple {U,n), hence, with the 
possibility of trajectories that are discontinuous at some point in the first component of 
the unknowns. 

Differentiating with respect to the equation defining the surface E and eliminating 
the variable m, we get the reduced two-dimensional system 

dU n-{go9)iU) dn 

Thanks to Q, function / is the second-order polynomial given in (15) and f'{U) = 
—2k{U — Uc)- The final step consists in introducing new coordinates to transform the 
system in the standard form ([s]). 

Introducing the variable (x, y) related with (f7, n) by 

U = U{x):=]^[[u]x + U.+U+] and n = n{y) ■.= ]^[[g]y + g{9.) + g{d+)] , 

system ( |17[ ) becomes 

^ ^ [a] G{x) - y dy ^ as k [u]^ _ 2n 
d^~ asK ^2 X ' d^~ 2[g] ^ ' 

where 

(18) G{x) := (2 (5 o o U)ix) - g{e.) - g{e+)) 

By rescaling the variable ^ by setting ^ = r we end up with the system ([s]) where 

(19) .:= ""-y . 

Both function G and parameter v depends on the values of the asymptotic states of the 
profile and thus can be regarded as functions of the parameters [u] and Uc- 
Let us also observe that the function 6 oU = {0 o U){x) is given by 

{eoU){x) = ]^{^-^Uc-u^u = ^{uc + \[u]x^ [u,-]^^[u]xy, 

thus, there holds 

-^(6* o U){x) = {x - Xc) where Xc := -- — > 0. 

ax 2K 7 [u\ 

As a consequence, if g is strictly increasing, the function G has the same monotonicity of 
9 o U and it is strictly increasing in [—1, 1] if and only if x > 1, that is if and only if 

I I 7 - 1 
\u+ — < Uc - 

7 

Otherwise, G is strictly increasing in [— l,Xc] and strictly decreasing in [xc, 1]- 
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Example 2.2. In the case g{0) = aO with o" > 0, the parameter v and the function G 
can be computed exphcitly 



Ta [9] 2(7-1)2 rcj 

2{0 o u){x) - 0+ - _ 

Thus, the structure is essentially the same of the Hamer model Q with a quadratic 



, 2{0 o u)ix) - 0+ - 7 2n 



function (7 (see Example 2.1). 



Higher powers in the temperature dependence of function g, i.e. g{9) = a 0^ with o" > 
and a a positive integer, give raise to complicate formulas, whose explicit expression is not 
particularly significant. It is relevant to observe that, being 0oU a homogeneous function 
of degree 2 with respect to [u] and Uc, the value [g] is homogeneous of degree 2a and thus 
v is homogeneous of degree 2(1 — a) and G is homogeneous of degree with respect to 
the same variables. In particular, when the ratio 

^ [u] ^1 U^-U+ 
Uc 2U- + U+ 

is kept fixed the function G does not change and the parameter v decreases as — [n] 
increases if a > 1. 

An eventual discontinuity in the trajectory keeps y fixed; moreover, if xn and Xr denote 



the values at the left and at the right of the jump, then and condition 06^ 

translates into 

(20) ^<-<»- 

Note that the value ({7+ + U- )/([/+ — U- ) is always strictly smaller than —1, so that jumps 
are always alloweded if x € [—1,0). 

Remark 2.3. The analysis is based on the choice ([2|, that playes a key- role in determining 



the specific expression (15) for the function /. For more general expression of pressure p 



and internal energye e, all of the subsequent analysis can still be performed if the function 



/, given by formula (14), is strictly concave. 

3. Analysis of the reduced system 

In this Section, we consider piecewise smooth solutions to the system 
dx _ 1 G{x) - y dy_l ^ ^ 
^^^> dC~ V X ' dC' 2^^^^ >■ 

with > and G smooth and such that G(ibl) = ±1. At any eventual jump point 1^, we 
assume that the couple {x^y) possesses left/right limits, denoted by (x^yi) and {xr,yr), 
and that 

xi < Xr, xe + Xr = 0, ye - Ur = 



System (21) has two singular points, P- := (-1,-1) and P+ := (+1,+1), and we are 



interested in the existence of a heteroclinic orbit connecting such that (x, y)(ibcx)) 
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P±. Since both critical points P± are saddles, such orbit is obtained by matching the 
unstable one-dimensional manifold of P_ with the stable one-dimensional manifold of Pj_ 



(see Lemma 3.1 ). 




Figure 1. The phase plane {x,y) with the graph of function G (continuous line) and 
the unstable/stable manifold of P- = (—1, —1)/P+ = (+1, +1) (dashed lines). 

The matching may happen either in a continuous or in a discontinuous way, the former 
case being related with the nature of the system at the singular line x = (see Lemma 



3.2) 



Lemma 3.1. For any z/ > 0, the critical points P± are saddles and, denoted by lA- the 
part of unstable manifold of P- in {x > —1} and by 5+ the part of the stable manifold 
of P+ in {x < 1}, there exist functions cp- = <^-(x) defined in [—1,0] and = 0+(x) 
defined in [0, 1] both differentiable and strictly monotone increasing, such that 

{ix,y)GU- : xG [-1, 0)} = {(x, y) : xG [-1, 0), y = ,^_(x)} 

{(x, y) G 5+ : X G (0, +1]} = {(x, y) : x G (0, +1], y = 0+(x)} 

Moreover, (j)+{x) < G{x) for any x G (—1,0) and 4>-{x) > G{x) for any x G (0, 1). 



Proof. The jacobian matrices of (21) at P± 

1 (±G'{±1) Tl 
u\ ^v'' 

have both determinant equal to —1; hence, both points are saddles. 

The manifolds U- and 5+ have tangent directions given by the vectors (/i" , v) and 
(— respectively, where 

„ 2u , _^ 2v 

^" G'(-l) + VG'(-l)2 + 4i/2' G'(+l) + VG'(+l)2 + 4z.2' 

Thus, if G'[—l) < 0, lies in the region {(x,y) : y > G(x)} in a neighborhood of P-; 
similarly, if G'{1) < 0, 5+ lies in the region {(x, y) : y < G(x)} in a neighborhood of 

If G'(±l) > 0, since fil < v/G'{-l) and < v/G'{+l), for x G (-1,1), the 

manifold lies in the region {(x,y) : y > G(x)} and the manifold 5+ in the region 
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{{x,y) : y < G{x)} in a neighborhood of P_ and P+, respectively. Since := {{x,y) : 
X E (— l,0),y > G'(x)} is positevely invariant and f]+ := {{x,y) : x G (0,1),?/ < G{x)} 
is negatively invariant, the U- is described in by a curve y = (p-ix) and the manifold 
5+ is described in $7-1- by a curve y = (j)-^-{x). Moreover, both functions are monotone 
increasing with respect to x since the corresponding derivatives, given by 

, , d , , , z^^ x(l — x^) 

(22) -Mx) 



dx^^' ' 2{G{(t>) 



do not vanish in the interior of Q- and $7+ respectively. Moreover, because of (22), the 
functions (p± cannot have vertical asymptotes in the regions (—1, 0) and (0, 1), respectively; 
therefore, they are well defined in (—1,0) and (0, 1), respectively. □ 

For a better understanding of the behavior of the curves at x = 0, we consider the 
half plane x > and rescale the independent variable by setting 

d d 

'^^^Tc = dv 



so that (21) takes the form 

(23) | = ^^('-^" 



Lemma 3.2. Assume G'(0) / 0,2^"^. The critical point Pq = (0,G(0)) of system ([23]) 
is an attractive spiral if G"(0) < 0, a repulsive spiral if < G'{0) < 2u^ and a source 
if 2v'^ < G'{0). In the latter case, there exists x^ G (0,1] and a monotone increasing 
function = ip-\-{x) defined in [0,X-|-] such that the fast unstable trajectory of Pq in the 
region = {x G (0, 1), y < G(x)} is given by {x,ip+{x)) . 

Proof The jacobian matrix at the singular point Pq := (0, G(0)) 

1 (G'{0) -1 
V \ v'^l2 

has real eigenvalues if and only if G'(0)^ — 2v^ > 0. In the opposite case, both eigenvalues 
have real parts of the same sign of G'(0). If 2v^ < G'(0)^, the eigendirections are given 
by (2/i^, z^) where 



± _ G^(0) ± ^/G'jO)^ - 2t/2 
^0 " 2^ 

Since 

< ^ < — < G'(0) 

the principal directions of the critical point (0,G(0)) both lies into the region y < G(x) 
for X > 0. The fast unstable direction, corresponding to the eigenvalue ^q, is uniquely 
determined and, while lying in Q^, is the graph of a monotone increasing function = 
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Remark 3.3. If the function G does not depend on the parameter v the vector field in 



(23) is a rotated vector field in the sense of Duff [3] in the regions 0±, as a consequence of 



the equaUty 



F G\ 1 , , 2^ 



(for the analysis of rotated vector field, see also [I^). Thus, the graphs of the function 
rotates anticlockwise as v increase and the map v ^ x'^ \s monotone decreasing with 
respect to v. Moreover, there exists a critical value vj^ such that x^ = 1 for v < vj^ and 
xK < \ ioi V > u^. 



From Examples 2.1 and 2.2, it is readily seen that the eventuality of a function G 
independent on v is rarely verified in concrete cases. Nevertheless, in the case of radiation 
hydrodynamics, the function G does not vary if the ratio 5 = —['u\/Uc is kept fixed. Thus 
the properties of the rotated vector fields return to be useful in the analysis of variations 
of the heteroclinic orbit connecting P_ and P+ when varying the asymptotic states U± 
keeping 5 fixed. 



A result analogous to Lemma 3.2 holds on the side x < 0, the main difference being 
related to the fact that the rescaling of the independent variable changes the orientation of 
the trajectories. In particular, we can state that, for 2z^^ < G'(0), there exists x_ G [—1, 0) 
and a monotone increasing function ij:- = tp- [x) defined in [x_ , 0] such that the fast stable 
trajectory of Pq in the region fi_ = {x G (0, 1), y < G{x)} is given by (x, 

Proposition 3.4. Let u > and G G C^([-l,l]) be such that G{±1) = 0. Then there 
exists a heteroclinic orbit connecting P_ with P^. Such orbit possesses at most a sin- 
gle jump, and such discontinuity is actually present if one of the following conditions is 
satisfied 

(24) G(0) < -1, or 1 < G(0), or G'{0) < 2u'^. 

Finally, ifG'{0) > the trajectory of the heteroclinic orbit is unique. 



Proof. 1. The trajectories {x,(p±{x)), defined in Lemma (3.1), converge to the y— axis at 
a finite value of C- Indeed, there holds 

° Jx{Co)dx/dC y^.(^Q) G(x) - </>±(x) 

The latter integral is finite if (j)±{0) / G(0). If 0±(O) = G(0), then (j)'^{0) ^ G'{0) as 
can be seen from the explicit expression of the principal directions of the singular point 



Pq of the rescaled system (23) (see Lemma 3.2). Thus, the function x/{G{x) — 4>±{x)) is 



integrable in small neighborhoods of 0. 

2. If (t>-{Q) = 4>+{^) = G(0), the conjunction of the two orbits 4>- and 0+ furnishes 
the heteroclinic orbit connecting P_ and P+. In this case, the profile of the corresponding 
traveling wave is continuous. 

3. To complete the proof of the first part of the statement, we assume that either 
(/)_(0) 7^ G(0) or 4>+{!^) 7^ G(0). Denoted by Tp±{y) the inverse functions of (j)± defined in 
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[— 1,0_(O)] and [(/)+(0),l], respectively, let us consider the function 

h{y) :=^+(y) + ^_(y) y G [<A+(0), <A_(0)]. 

Then, by the properties of (j)±, the function h is strictly increasing and 

h{MQ)) := i^-iMO)) < < ^+(0-(O)) = h{(t>-m. 

Hence, by continuity of h, there exists a single value yc such that h{yc) = 0. The connection 
of the trajectory cp^, truncated at x = ip-{yc), and 0+, truncated at x = il^j^iyc) gives the 
desired heteroclinic orbit. 



4. Any of the conditions in (24) implies that one of the trajectories 4)± does not pass 
through the point (0,G(0)). Indeed, if G(0) < —1, since the function is monotone 
increasing, then (/'-(O) > — 1 > (5(0). Similarly, if 1 < G(0), then i?I)+(0) is strictly greater 



than G(0). Finally, if G'(0) < 2z^^, the point Pq is a spiral for the rescaled system (23) 
and no monotone trajectory is alloweded to get onto it. 

5. To complete the proof, it is sufficient to observe that more can one jump may appear 



only if there are smooth solutions of (21) passing from the positive side {x > 0} to the 
negative side {x < 0}. Because of the structure of the vector field defining the system, the 
corresponding eventual trajectories should pass through the point Pq = (0,G(0)), but, in 
the case G'(0) > 0, the singular point Pq is repulsive on {x > 0} and attractive on {x < 0} 
so that no transition from positive to negative values of x is possible. □ 



Remark 3.5. In the case G'{0) < 0, the system (21 ) possesses trajectories passing through 
the point (0,G(0)) and going from the region x > to the region x < 0. Thus, in 
principle, solutions may have more than a single jump point and we cannot exclude that 
there exists also other possible orbits with the same asymptotic states. In the concrete 



case of application of Proposition 3.4 to ^ and ([T]) the condition G'(0) > is satisfied. 



4. Back to radiation hydrodynamics 



Next, let us analyze the system (21 ) and the properties of the corresponding heteroclinic 
orbits in the case of radiation hydrodynamics as derived in Section [2j For the reader 
convenience, let us recall the formulas 

G(x) ■.= J-{2{goeo U){x) - g{e^) - g{e+)) and - - t^]' 



[^p-.. ^^^j 



(see formulas (18)-(19)). As stated in the Introduction, we introduce the key-parameter 



5 := — [u]/f7c) so that the expression of the compound function (OoU) = {9oU){x) becomes 

(25) ^0oU){x) = ^[l-\5x^ {l + \l6x 

that describes the values of temperature as a function of the auxiliary variable x. 
We concentrate on a series of significant issues: 

- to determine threshold for the presence/absence of disconinuities; 

- to analyze the monotonicity of the temperature profiles; 
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- to study the regimes 5 — )• and J — )• 2/7; 

- to show a short coUection of numerical experiments. 

At the end of the Section, for completeness, we also mention how the properties of the 



reduced system (21) translate in the case of the Hamer model. 



Discontinuous profiles. As stated in Proposition |3.4[ profiles of the traveling waves are 



discontinuous if condition (24) is satisfied. If g is monotone increasing and [g] > 0, there 



hold G(0) < -1 if and only if {6 o U){0) < 6- and G(0) > 1 if and only if 0+ < (0 o U){0). 



Taking advantage of the expression (25), the previous relations can be rewritten as 



G(0)<-1 ^ l<{l + l6)(l-l^5] ^ 6<^^'-^^ 



2 J \ 2 ' J - 7 

G(0)>1 « (l4*Vl + ^*)<l ^ *>'*^-'' 



2 J \ 2 ' J - - 7 

For 7 > 1 and > 0, the first condition is never satisfied. 
Additionally, there holds 



7-1 [u] dg ( U[ 



so that condition G'(0) < 2z^^ is equivalent to 

, , Mo,dg(U^\ [uf , 2(7-1)3 

(26) -X^Uci^i^] < - Vr where Mq — ^' ' 



A^^de\-iRj [g] i?7(7 + l)2cT2 



For g{9) = aO, condition (26) becomes 



. mi ^ \/2(7 - 1)2 err 

> — - where mi := — — — 

A ^7(7 + 1) 

For g as in ^ with a integer greater than 1, the explicit formulas corresponding to the 

condition ( |26[ ) become very intricated and they do not seem to be significant. In any 

case, let us stress that differently with respect to the condition G(0) > 1, that depends 

only on the relation between the parameter 6 and the value of the constant 7, constraint 

G'{0) < 2i/2 depends also on the other characteristic parameters of the system. 

Summarizing, we can state that the heteroclinic orbit is discontinuous if 

(27) 6 > (5ju^p := 2^1^ 

7 

Note that, if 7 > 2, the region in the plane (C/c, Kti]]) for which this inequality is satisfied 
does not intersect the attainable region 5 < 2/7. 



Temperature spikes. The temperature profile is given by (25) calculated at x = x(C), 
the first component of the heteroclinic orbit of the reduced system (21). Since x = x{Q 
is strictly increasing, the monotonicity of the temperature is controlled by the sign of its 
first derivative 

^{0 oU){x) = -^^{x - Xc) where Xc:=- — ^\ 
dx 2R 70 
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Hence, the profile has a change in monotonicity if and only if 

7 - 1 

(28) 6 > (5spikc := 

7 

Such absolute maximum for the temperature is sometime referred to as the ZeI'dovich spike. 



Note that the condition (28) is optimal and depends only on the value of the constant 7 



(in particular, it does not depend on the specific form of the function g\). For 7 > 3, the 



condition (28 ) is not satisfied if <5 < 2/7, meaning that for such values of 7 the temperature 
profile is always monotone increasing. 

The location of the absolute maximum point for can be either at the jump point, if 
present, or in a region of regularity of the profile. In the latter case, the maximum of the 
temperature is attained at x = Xc and its value is 

(29) ^_:=(0ot/)(x,) = ^^^^ 

If the point of change of monotonicity coincides with the jump point, the values ^max gives 
only an estimate from above of the maximum values for the temperature. It is particularly 
relevant to note that the expression of 6'„iax does depend only on 7, R and Uc- 

Estimate ( |29[ ) can be compared with the corresponding one determined in [13] (formula 
(104.67), p.573) for 7 < 3 

:= (3-7) 0+. 

By using the expression for 9^ in terms of Uc and [u], we obtain 

gmax (7+1)' 1 

9m 47(3-7) (1-<5/2)(1 + 7^2)' 
Since 5 £ (0,2/7), we infer 

(30) _J_^0^ (7 + 1)^ 

^ ' 3-7 - - 8(3-7)(7-l)' 

In the significant cases 7 = 5/3 and 7 = 7/5, the above estimates becomes 

5 3 ^max ^ -I ^ ^ ^ ^max ^ 9 

T = — : — < < 1 'y = — : — < < — . 

'3 4-^M~ 5 8- 9m ~ 8 



The expression (30) can be rewritten as 



1_1^<^<1+ (37-5)^ 



3-7 " - 8(3-7)(7-l)' 
thus, in the regime 7 < 2, the first term in this equality chains is strictly smaller than 
1 and the last term, apart for the (physical!) case 7 = 5/3, the estimate from above is 
strictly greater than 1. 



Limiting regimes. The expression (25) shows that 9 oU has the form C p{x;j,6) with 
C = Uc/'J R and p a polynomial of degree 2 in x with coefficients depending only on 6 
and 7. Therefore, in the case of functions g with the power- law form described in ([2]), the 
function G in the reduced system ^ does not changes if the parameters 7 and 6 are kept 
fixed. Here, we want to analyze the behavior of the heteroclinic orbits built in Section [3] 
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for 5 prescribed in the regimes 5 — )• and — )• 2/7. To this aim, we need to investigate 
the function G for such values of 5 and analyze the variations for the parameter v. 
As 5 — )• 0, formula (|25|) can be rewritten as 



[9oU){x) = ^{l + \{^-l)5x^+o{5). 

and substituting in the definition of G, we infer 

G{x) = x + o{l) 5^Q. 

Thus, the function G is strictly increasing and the temperature profile is monotone in- 
creasing. As 5 — )■ 2/7, the relation (25) becomes 

{9 o U){x) = (7 - x) (1 + x) + o{5 - 2/7) 
7^ H, 



Thus, we have 



where 



= G^,„(x) + o(,5-2/7) 
G^,«(x):=-l + 2 ^^""^ ^^"^ 



2 

1 

7 



,7-1 2 

In this case, the function G changes its monotonicity if and only if 1 < 7 < 3 and, if this 
is the case, has an absolute maximum point at (7 — l)/2; moreover, G(0) > 1 if and only 
if 7 < 2. In particular, the temperature profile exhibits the Zel'dovich spike. 

Next, we examine the parameter u and its behavior for 5 fixed and either [u] — )• or 
[u] — )• —00 (equivalently, either C/c — >• or C/c — >• +00). The assumption g{9) = a 9^ carries 
also the identity 

[g] = a{9+ - 9^){9X-^ + ■■■ + 9^-^) = -Ga[u] f/^""! 

for some positive constant C only on 7, 5, R and a. Inserting in the definition of given 
in (19), we end up with 

_ C5 _ C(^2a-l 

for some positive constant C dependent only on 7, A = p± U±, a, aa, (Jg- Hence, for a > 1, 
V — )• +00 as [u] — )• (or Uc — )• 0) and 1/ — )• as [u] — )• —00 (or Uc — ?• +00). The case a = 1 
is different, since v turns to be constant for 5 fixed (see Example 2.2). For this reason, we 
only consider the case a > 1. 

As [u] —7- +CXD, in an appropriate scale, system ([s]) formally reduces to 

dx ^ G{x) - y ^ = Q 
dC, X ^ dC, 

For (5 —7- 0, such structure forces the unstable manifold lA^ of P_ = (—1, —1) and the stable 
manifold 5+ of P_ = (+1, +1) coincide with the graph of the function G{x) = x + o{5). 
The corresponding temperature profile is continuous and strictly monotone increasing. 

For S — )• 2/7, and 7 < 2, then the unstable manifold U- of P_ = (—1, —1) coincides 
with the graph of the function G-y^a] the stable manifold 5+ of P- = (+1, +1) is given by 
a horizontal segment at height y = 1. The heteroclinic orbit is obtained by the connection 
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of such manifold at the (unique) point where the stable manifold U- is at height 1 (see 
Figure [2| left). Since G^,a is symmetric with respect to Xc = 5(7 — 1), the value x_ of 
the jump is equal to 2xc — 1 = —(2 — 7) and the corresponding value is 2 — 7. The 
maximum of the temperature profile is attained at the jump point if and only if x+ is 
greater than or equal to Xc, that is if and only if 7 < 5/3. 




■1 -0.5 as 1 -1 .0.5 as 1 




_j I I i l_i , , I I I , , ■ I ■ t ■ : I ; i—i I I ■ ■ i ■ I L 1 ' I I I I 1 I 1 I I I 1 I I I I I I I I I I I , I l_ 

10000 £0000 30000 40000 50000 60000 20000 40000 60000 60000 100000 120000 140000 



Figure 2. Re gime 5 — >■ 2/7 and [it] — >■ — oo. Upper line: the phase plane {x,y) with 
the graph of the function G (continuous line) and the heteroclinic orbit connecting the 
points P_ = (—1,-1) and P+ = (+1,+1) (dashed lines for the continuous part and 
dashed-dotted line for the jump). Case 7 < 2 (left) and case 7 > 2 (right). Lower line: 
corresponding temperature profiles for the case 7 < 2 (left) and case 7 > 2 (right), both 
with an inside box showing a zoomed version of the spike. 

For 5 — )• 2/7, and 7 > 2, then the unstable manifold U- of P_ = (—1, —1) coincides 
with the graph of the function G^^a on the interval [—1,0]. The stable manifold 5+ of 
P_ = (+1,+1) is given by the union of the horizontal segment connecting the point 
(7 — 2, 1) and (1, 1) and the graph of the function G^y^^ for a; G [0,7 — 2] (see Figure [2| 
right). 

For [u] — )■ 0, again in an appropriate scale, system (|3]) becomes 

In this regime, the limiting dynamics is very simple: the unstable manifold U- of P = 
(— 1, —1) coincides with the vertical half-line {x = —l,y > —1} and the stable manifold 
5+ of P_ = (+1,+1) with the vertical half-line {x = l,y < 1}. The heteroclinic orbit 
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is described by a direct jump from U- to 5+. The corresponding temperature profile is 
monotone for 5 — )• 0^ and non-monotone for 5 — )• 2/7. 

The intermediate cases for u can be deduced quahtatively by noting that, for 5 fixed 
and V varying, the system ([s]) is a rotated vector field, as already observed in Remark 3.3 



In particular, as v increases, that is as \ decreases, the trajectory determining the 

manifolds lA^ and 5+ rotate counter clockwise passing from one limiting configuration to 
the other. 

Numerical experiments. One of the advantage of the reduced system ([3| resides in 
its semplicity and in the possibility of a numerical approximation of the structure of 
the radiative profiles by means of a standard solver for ordinary differential equations, 
complemented with the conditions relative to the eventual jump point. We present here a 
series of experiments where parameters have been chosen only in part realistic. The value 
for the adiabatic constant 7 is taken equal to 5/3 and the constant R is set equal to 8.31. 
All of the constants relative to the radiative coupling, a, cts and r, are chosen equal to 1. 
The function g has the form expressed in Q. As quoted, the realistic exponent a is 4; 
nevertheless, for computational advantage, we choose to set a = 2 (the case a = 1 has 
different behavior in the limiting regimes). Also, to reduce the number of free parameters, 
we fix p_ at the value 0.1. 

Finally, the different cases depend on the choices of the ratio 5 and values of Uc (the 
average of the velocities U± at the left and at the right of the profile). Since 7 = 5/3, the 
thresholds for 5 are 

2 

(^spike = 0.4 (^■jump = 0.8 - = 1.2 

Being interested in the behavior for weak/strong radiation and weak/strong shocks we 
choose 5 = 0.6 (Figure |4| left column) and (5 = 1.0 (Figure |4| right column) and Uc equal 
to 10,50 and 100. 

Coeherently with the result previously deduced, all of the profiles are non-monotone 
since both the values of 5 are taken above the threshold dspike = 0.4. The choice 5 = 1.0 
is also above the threshold Jjump = 0.8 and thus the profile are discontinuous. For all of 
them, the spike coincide with the jump point. The choice 6 = 0.6 is consistent with both 
continuous and discontinuous profiles and this is revelad by the fact that for Uc = 10, the 
temperature has a jump, and for both Uc = 50 and Uc = 100 the profile is continuous. 

Higher values for 6, closer to the limiting value 2/7, give raise to pictures resembling 
the one already presented in Figure [2] (left column). In particular, the width of the spike 
region shrinks to a single point. 

Hamer model. In the case of the system Q with f(s) = ^ s^, if g is monotone increasing, 
G{0) < — 1 if and only if c > n_ and G(0) > 1 if and only if < c. Since n-(_ < c < n_ 
for any admissible shock wave, discontinuous orbits may appear only as a consequence of 
the condition G'{0) < 2z^^, that, in the present case, translates into 

du ^2 [g]- 
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Figure 3. Temperature profiles in tiie cases S — 0.6 (left column) and S = 1.0 (right 
column), corresponding to the values Uc ~ 10, 50, 100 (increasing from top to bottom). 



Considering g{u) = a u, a > 0, (see Example ( |2.1[ )), the above condition becomes 

|m+ — M- I > 

that describes the (sharp) threshold for discontinuous profiles (see [8]). 
The case g{u) = av? gives the condition 

a c 

showing that discontinuous profiles are possible also in the regime of small — u_ | , if the 
sonic value c is sufficiently small. Note that in this case, the function G is non-monotone 
in [—1, 1] if and only if 

|m+ — M_| > 2 |c|, 

hence, if u > l/\/2, both non-monotonicity of G and continuity of the profile are compat- 
ible for some appropriate choice of the asymptotic states u± . 
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5. Conclusions 

The analysis performed shows that shock waves are robust patterns for the system ([T]) 
under the assumption Q. Indeed, the existence of such structures is not hmited to specific 
ranges for the asymptotic states; oppositely, they do exist for any regime consistent with 
the reduced hyperbolic system, obtained by disregarding the radiation effects. Addition- 
ally, the system of algebraic-differential system can be dealth with in a rigorous manner 
with a restricted number of technicalities and, at the same time, with a good number of 
quantitative informations relative to the internal structure of the shock transition, specif- 
ically with respect to the presence of change of monotonicity. This allows to extend the 
analysis carried out in the region 5 —t- 2/7 in [13] to any regime and gives a sound basis 
to the result formally derived in [7|. Also, the observation that, for fixed 6, the reduced 
system ([s]) depends on the value \U^ — U-\ (or, equivalently, on Uc) as a rotated vector 
field permits, in principle, to follow the changes in the qualitative properties of the profiles 
in a "monotone" way with respect to the variation of such value. 

The relative simplicity and solidity of the structure of the problem is also revealed from 
the fact that there is a great freedom in the possible form for the coupling term g -though 
the power-like form is the most accredited version- and that there is a reduced number 
of parameters required to classify different forms for the profiles. More general forms for 
pressure and internal energy should also preserve the existence of the radiative profiles so 
that it is reasonable that the results presented in this article extend to more general fluids. 

A natural subsequent step is the stability analysis for such radiation profiles, in partic- 
ular in the presence of spikes. Such problem could be approached taking inspiration from 
other classes of systems of partial differential equations in the context of conservation laws 
that support non- monotone traveling wave as in the case of combustion waves. Also, it 
would be very interesting to analyze the case of radiation hydrodynamics in more general 
forms and specifically without taking advantage of the assumptions that reduce the model 
to a hyperbolic-elliptic system of partial differential equations. In particular, considering 
the native form of a continuous-kinetic coupled model would be a very challenging and 
stimulating research project. 
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